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ABSTRACT 

We present the hydrodynamic BL Herculis-type models which display a long-term 
modulation of pulsation amplitudes and phases. The modulation is either strictly pe- 
riodic or it is quasi-periodic, with the modulation period and modulation pattern 
varying from one cycle to the other. Such behaviour has not been observed in any 
BL Her variable so far, however, it is a common property of their lower luminosity 
siblings - RR Lyrae variables showing the Blazhko effect. These models provide a 
support for the recent mechanism proposed by Buchler & Kollath to explain this still 
mysterious phenomenon. In their model, a half-integer resonance that causes the pe- 
riod doubling effect, discovered recently in the Blazhko RR Lyrae stars, is responsible 
for the modulation of the pulsation as well. Although our models are more luminous 
than is appropriate for RR Lyrae stars, they clearly demonstrate, through direct hy- 
drodynamic computation, that the mechanism can indeed be operational. 

Of great importance are models which show quasi-periodic modulation - a phe- 
nomenon observed in Blazhko RR Lyrae stars. Our models coupled with the analysis 
of the amplitude equations show that such behaviour may be caused by the dynamical 
evolution occurring in the close proximity of the unstable single periodic saddle point. 

Key words: hydrodynamics - methods: numerical - stars: oscillations - stars: vari- 
ables: BL Herculis - stars: variables: RR Lyrae 
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1 INTRODUCTION 

Nonlinear modelling is one of the key methods to study the 
large amplitude variability of classical pulsators: RR Lyrae 
stars and Cepheids. One of the notable successes of non- 
linear pulsation theory is the explanation of resonant ef- 
fects in classical pulsators. One of the effects is the bump 
progression in the light/radial velocity curves of clas- 
sical Cepheids pulsating in the fundamental mode (so- 
called Hertzsprung progression). The effect is caused by 
the 2:1 resonance between the fundamental mode and 
the second overtone, the latter mode being resonantly 



excited to high amplitude (e.g. Simon fc Schm idt 



1976 



1990 



Kovacs fc Buchlerlll989l: iBuchler. Moskalik fc Kovacsl 
Bono. Marconi fc Stellingwerf 2000). Due to resonant non 



linear phase synchronisation the second overtone is not vis- 
ible separately, but manifests itself as a distortion in the 
light/radial velocity curve. A similar effect is present in 
classical Cepheids pulsating in the first overtone (see e.; 
iKienzle et all Il999l : iFeuchtinger. Buchler fc Kollathl [2000 
Another interesting resonant effect, which we understand 
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thanks to nonlinear pulsation theory, is a period-doubling 
behaviour - alternating deep and shallow minima in the 
light/radial velocity curves. It is a characteristic feature of 
RV Tau variables, a group of pulsators o ften considered to 
be a subgroup of type-II Cepheids (e.g. IWallersteinl 120021 : 
ISzaba dos 2010). The period doubling e ffect is caused by the 
half-in teger resonances, as analysed by Moskalik fc BuchlerJ 
(1990). In case of RV Tau variables, the 5:2 resonance be- 
tween the fundamental mode and the second overtone is 
crucial. Period doubling was also discovered in Blazhko 
RR Lyrae stars observed with the satell ite mission Kepler 
dKolenberg et al.ll2010l: ISzabo et alj|2010h and explained by 
iKollath. Molnar fc Szabd i|201ll ) as a result of a 9:2 reso- 
nance between the fundamenta l mode and the 9th over- 
tone. Recently ISoszvhski et al.l |201ll ) and ISmolec et al.l 
(2012) reported the discovery of the period doubling effect 
in a BL Her star, demonstrating the predictive power of 
the nonlinear pulsation theory. Existence of period-doubled 
BL Her stars was in f act p redicted t wenty years earlier 
by IBuchler fc Moskalikl l|l992T l (see also IBuchler fc Buchlerl 
1 1994 ). who found the e ffect i n their radiative hydrodynamic 
models. ISmolec et all (|2012h confirmed that the 3:2 reso- 
nance between the fundamental mode and the first over- 
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tone is the cause of t he alt ernations, as analysed earlier by 
iBuchler fc Moskalik! (fl992t ). 

In lSmolec et alj (I2012T ) we have studied the period dou- 
bling effect in BL H er models with our state-of -the-art con- 
vective hydrocodes |Smolec fc Mos kalik 2008a). During the 
test computations with the decreased eddy viscosity we iden- 
tified a class of models showing the modulation of pulsa- 
tion amplitude and phase, which is either strictly periodic 
or quasi-periodic. Such behaviour has not been detected in 
any BL Her star so far. Also, because of the strongly reduced 
eddy- viscosity, the pulsations are violent and spurious spikes 
appear in the light curves. Still, the models are of great im- 
portance for the lower luminosity cousins of BL Her stars 
- RR Lyrae pulsators, in which the more or less periodic 
modulation of pulsation amplitude and phase - the Blazhko 
effect - is a common property (e.g. IKolenbe rg 2008). 

The amplitude modulation in RR Lyrae stars is one 
of the most disturbing problems of stellar astroph ysics. Al- 
thoug h it was discovered more than century ago (|Blazhkd 
1907), its origin is still mysterious. During the recent 
years extensive ground-based observation campaigns (e.g. 
IKolenberg et al.l 120061 ; lJurcsik et alj l2009h and a top qual- 
ity and nearly continuous observations of space missions 
CoRoT and Kepler allowed to rule out the two models pro- 
posed to explain the Blazhko effect, the magnetic oblique 
rotator/pulsator model and the non-r adial resonant rota- 
tor/pulsator model f see iKovacsl [2009 . for a review). One 
of the main drawbacks of these models is the predicted 
clock-work modulation, while it became clear, with the ad- 
vent of satellite data, that the B lazhko effect can be a very 
irregular phenomenon (see e.g. iGuggenberger et"afl l201l|. 
120121; IKolenberg et all I20T1I ). The recent idea of IStothera 



20061 ). which ssumes that modulation of turbulent con- 
vection by transient magnetic fields causes the modula- 
tion of pulsation was also que stioned (|Smolec et al.l l201ll ; 
iMolnar. Kollath fc Szabdl I2012T ). The discovery of period 
doubling effect in some Blazhko variables led to the new 
model behind the Blazhko mo dulation - the radi al mod e 
resonance model proposed by IBuchler fc Kollatbl (|201ll ). 
Using the amplitude equations' formalism (AE, see e.g. 
IBuchler fc Goupillll984l l IBuchler fc Kollatbl l|201ll ) showed 
that the same resonance that causes the period doubling in 
Blazhko variables, i.e. the 9:2 resonance between the funda- 
mental mode and the ninth overtone, can also cause either 
periodic or chaotic modulation of the pulsation. The AE for- 
malism is a powerful tool to study the nonlinear pulsation, 
specifically the possible limit cycles and their stability. The 
nature of the computed limit cycles, however, depends on 
the arbitrarily assumed values of the saturation and reso- 
nant coupling coefficients, f or which realistic values are too 
difficult to compu te (see e.g. lKlapp. Goupil fc Buchlerll 19851 : 
iNowak owski 2005|) . Consequently, it is not clear whether the 
specific solution of the AEs in which pulsation is modulated, 
can be represented in the real stars. This must be confirmed 
with realistic nonlinear hydrodynamic models. Attempts to 
find modulation of pulsation in hydrodynamic models of 
RR Lyrae stars did not lead to success so far, likely be- 
cause of the surface nature of th e involved 9th order over - 
tone, which is difficult to model (|Buchler fc KollathllioTl} ). 
In this paper we present the hydrodynamic models of BL Her 
variables, more luminous, but otherwise (masses, chemical 
composition, evolution history) siblings of RR Lyrae stars. 
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Figure 1. The HR diagram for the computed models of M = 
0.55 Mq and Z = 0.0001. Instability strips for the fundamental 
mode and the first overtone are indicated with the light- and dark- 
shaded areas respectively. Loci of several resonances are plotted 
as indicated in the key. Horizontal line segments show the domain 
of modulation of the pulsations. 



In these models a low order 3:2 resonance, between the fun- 
damental mode and the first overtone causes the modula- 
tion of pulsation. Our calculations demo nstrate for the first 
time that the mechanism proposed by IBuchler fc Kollath] 
l|201ll ) can be indeed operational in hydrodynamical mod- 
els. In addition, they illustrate how the quasi-periodic mod- 
ulation may arise. Thus, they provide a strong support for 
the radial mode resonance model of the Blazhko effect. 

We first present our hydrodynamic models which dis- 
play modulation of pulsation (Section [2]) and next we 
demonstrate that the computed behaviour can be captured 
and understood with the AE formalism (Section [3J. In Sec- 
tion[4]we comment on the reliability of the models and mech- 
anisms responsible for the wealth of detected behaviours. 
Finally, in Section [S] we discuss the implications for under- 
standing the Blazhko effect in RR Lyrae stars. 



2 HYDRODYNAMIC MODELS 

The hydrodynamic models were comput ed with the Warsaw 
nonlin ear, convective, pulsation codes l|Smolec fc Moskalikl 
2008a). Numerical parameters (zoning) and convective pa- 
rameters, ex c ept o f eddy viscosity, are the same as in 
I Smolec et al.l (|2012|) (see their section 3.1 and table 2, set 
PI). The eddy viscosity (a m parameter) was decreased from 
0.25 to 0.05. The resulting fundamental mode instability 
strip (IS) - see Fig. [T] — is very wide, the red edge was not 
reached in our computations. For the first overtone we ob- 
serve an additional thin instability strip extending to the 
red of the usual first overtone IS and towards the higher 
luminosities. Consequently, our models do not present a re- 
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Figure 2. Radius variation for the three models of 136 Lq and 
different effective temperatures. 



alistic scenario for BL Her pulsation (see Section [4] for de- 
tailed discussion). Yet, they show a very interesting feature, 
namely the periodic amplitude modulation on top of the pe- 
riod doubling alternations - see bottom panel of Fig. [2] for 
an example. In Fig. [J] we show the location of these inter- 
esting models in the HR diagram (horizontal line segments). 
Over-plotted are loci of several half-integer resonances. Each 
of them can, in principle, c ause the period doubling e ffect 
and modulation of pulsation (|Moskalik fc Buchlerll990l . Sec- 
tion[3]). Although the modulation domain coincides with the 
loci of the 7:2 resonance with the 4th overtone, the modula- 
tion of pulsation is most likely caused by the 3:2 resonance 
with the first overtone, which is r esponsible for the p eriod 
doubling effect in BL Her models (|Smolec et al.ll2012T ). The 
origin of the modulation will be discussed in more detail in 
Section 2] 

Because of strongly decreased eddy viscosity the models 
are numerically demanding. Violent pulsation require very 
short time steps leading to extremely slow integration. Con- 
sequently, only a limited nonlinear model survey was con- 
ducted. All models have M = 0.55 M and Z = 0.0001 (X = 
0.76). We used t he OP opacities ( SeatonlfeoOSl ') and adopted 
a solar mixture of lAsplund et al.l ( 2004 ) . The nonlinear mod- 
els were computed along six lines of constant luminosity, 
L e {120L Q , 128L Q , 136L , 144L , 152L Q , 160L Q } 
and cover a limited range in effective temperatures - large 
enough to determine the boundaries of the domain in which 
modulation of pulsation takes place. We first describe the 
model sequences in which all models show strictly periodic 
modulation of pulsation (Section l2.ll) . Next we describe the 
model sequences along which models with quasi-periodic 
modulation were found (Section 12. 2[) . 



2.1 Models with strictly periodic modulation of 
pulsation 

Strictly periodic modulation was found in all models along 
four sequences of highest luminosity (see Fig. [T]). The pul- 
sation scenario is discussed based on the model sequence of 
L — 136 L©. The radius variation for the three models of 
different effective temperatures is presented in Fig. [2] The 
hottest model (T e g = 6300 K, top panel of Fig. ^ displays 
a single-periodic, fundamental mode pulsation. As we de- 
crease the model's effective temperature the period-doubling 
domain emerges. The radius variation is illustrated in the 
middle panel of Fig. [2] for model of T cft = 6180 K. Alter- 
nation of radii maxima and minima is clear. At still lower 
effective temperatures, periodic modulation of pulsation on 
top of the period doubling effect is found in the models - 
bottom panel of Fig. [2] The domain in which we detect this 
behaviour is not large, its width is « 50 K (see Fig.[T]). To the 
red of this domain the period doubling only (no modulation) 
behaviour is present again, and next we find a single-periodic 
pulsation. 

The described behaviour is also clear from the analysis 
of the frequency spectra of the three models of Fig. [5] pre- 
sented in Fig. [3] For better readability the frequencies are 
normalised with the fundamental mode frequency, /o. In the 
top two panels of Fig. [3] we show the frequency spectra for 
the two hotter models for which no modulation was detected. 
For single periodic pulsation in the model of T e g — 6300 K 
(top panel of Fig. [3} only the fundamental mode and its 
harmonics are present in the frequency spectrum. For the 
cooler model of T G g = 6180 K (middle panel of Fig. strong 
signal at the half-integer frequencies (sub-harmonics) is de- 
tected as well. It is a characteristic signature of a period 
doubling effect. After prewhitening with the fundamental 
mode and its harmonics (T e g = 6300 K) and sub-harmonics 
(T c ff = 6180 K) we find no additional signal in the data. 
In the bottom panel of Fig. UJ we show the prewhitening 
sequence for the model of T c g = 5892 K, that shows the 
modulation of pulsation on top of period doubling (see bot- 
tom panel of Fig. 01 ■ The top panel in the prewhitening se- 
quence shows the frequency spectrum of the original hydro- 
dynamic time-series. The fundamental mode, its harmonic 
and two sub-harmonics are well visible. At /o, l/2/o and at 
low frequencies, additional peaks are also present, which be- 
come clearly visible after prewhitening with the fundamen- 
tal mode its harmonics and subharmonics (middle panel in 
the prewhitening sequence). These peaks correspond to the 
modulation of pulsation. With the removed frequencies they 
form the equidistant multiplets. The separation between the 
multiplet components corresponds to the modulation fre- 
quency, / m , which is also detected independently in the fre- 
quency spectrum. For the discussed case the modulation pe- 
riod is P m = l//m = 63.8 days, which corresponds to 27.6 
fundamental mode pulsation cycles. 

The modulation of pulsation can also be illustrated 
through construction of the return maps. In Fig. [4] we show 
the first return map for the radii maxima, i.e., a plot of 
-Rmax vs - Chi where n counts the pulsation cycles, for 
the three discussed models. The return maps are plotted 
for the full model integration, including the initial transient 
phase. For single periodic model of T e g = 6300 K, after ini- 
tial transient, the consecutive radii maxima fall at a single 
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Figure 3. Frequency spectra for the models of Fig. [2] Top panel: 
frequency spectrum for model of T B g = 6300 K that exhibits single 
periodic pulsation; middle panel: frequency spectrum for model of 
T B ff = 6180 K that exhibits period doubling effect; bottom panel: 
Prewhitening sequence for model of T e g = 5892 K with periodic 
modulation of pulsation. Removed frequencies are marked with 
vertical dashed lines 



point located on the diagonal (cross in the leftmost panel of 
Fig. [jQ. It corresponds to a fundamental mode limit cycle 
pulsation. For models of T cff = 6180 K and T cff = 5892 K, 
which show the period doubling effect, it is convenient to 
connect each second pair of radii maxima (middle and right 
panel of Fig. After initial growth of the amplitude of 
single-periodic pulsation the so constructed curves diverge 
as a result of the raise of the period doubling effect. At the 
limit cycle pulsation the return maps for the two models are 
different. For the model of T c g = 6180 K, that shows the 
period doubling effect and no modulation (middle panel of 
Fig- HI) the consecutive radii maxima alternate between two 
points, located symmetrically with respect to the diagonal. 
For the model of T c s — 5892 K that shows the periodic mod- 
ulation of pulsation (right panel of Fig.|3}, each second pair 
of radii maxima follow a route along a separate loop. The 
modulation is strictly periodic, each loop can be described 
by a single curve - the broadening visible in Fig. U is caused 
by connecting with straight line segments. 

The described properties of modulation are qualita- 
tively the same for model sequences of higher luminosity, 
i.e. for L G {144 Lq, 152 Lq, 160 Lq}. In the upper panel 
of Fig. [5] we plot the modulation periods as a function of 
model parameters: luminosity and effective temperature. In 
the bottom panel the corresponding modulation amplitude, 
i.e. the amplitude of the highest modulation peak at the 
vicinity of fo (either /o — f m or /o + / m ) is plotted. The typ- 
ical modulation periods vary from ~ 16 to ~ 35 fundamental 
mode pulsation cycles (45 - 106 days) . In general, the higher 
the model luminosity and the higher its effective tempera- 
ture the shorter the modulation period. The modulation am- 
plitudes decrease towards the edges of the domain in which 
modulation of pulsation was found. This indicates that the 
appearance of the modulation domain, as the model's effec- 
tive temperature is changed, is a supercritical bifurcation. 
Also, for L ^ 136 Lq, the higher the luminosity the lower 
the modulation amplitudes. It is clear that the modulation 
domain shrinks towards higher luminosities and likely does 
not extend far beyond 160 Lq . 

In case of the least luminous models of 120 Lq and 
128 L© sharp features are clearly visible in the progressions 
of modulation periods and amplitudes with effective tern- 
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Figure 5. Modulation periods in units of the fundamental mode 
pulsation period (upper panel) and modulation amplitudes (bot- 
tom panel) as a function of effective temperature and luminosity 
(different line types). 



perature. These models are discussed in detail in the next 
Section. 



2.2 Models with quasi-periodic modulation of 
pulsation 

In Fig. [6] we plot the return maps for the limit cycle pulsa- 
tion (initial transient phase was omitted) for several mod- 
els along the 128 Lq sequence. As before, each second pair 
of points were connected (period doubling). For the coolest 
and two hottest models the return maps are qualitatively 
the same as discussed in the previous Section (two separate 
loops). These models show strictly periodic modulation of 
pulsation. To explain the sharp increase of the modulation 
periods visible in Fig. [S] we focus our attention on the three 
models with effective temperatures of 5860 K, 5872 K and 
5876 K. Their modulation periods increase and are equal to 
68.5 d, 73.6 d and 106.7 d, respectively. The radius variation 
for these models is plotted in Fig. [7] The return maps have 
a different topology from that discussed so far. They have a 
shape of figure '8'. For the two models of 5860 K and 5872 K 
the pulsation modulation is periodic, but the higher/lower 
maximum pulsation cycles trade places from one modula- 
tion cycle to the other. This behaviour can be followed in 
Fig. [7] (in particular in rightmost panels). Consequently, the 
pattern of period doubling repeats after two modulation cy- 
cles. We note that this behaviour was recently found in the 
Kepler short cadence data for RR Lyrae (Molnar et al., sub- 
mitted to MNRAS), although the swaping of higher/lower 
maximum cycles is very irregular. The frequency spectrum 
for the model of T e s = 5872 K (it is qualitatively the same 



for the model of T e g = 5860 K), after prewhitening with the 
fundamental mode and its harmonics, is plotted in the top 
panel of Fig. [5] At the position of the sub-harmonics, i.e., 
exactly at l/2/o and 3/2/o, we detect no signal. Instead a 
comb of equidistant frequencies, with peak separation cor- 
responding to / m is centred exactly on the location of the 
missing sub-harmonics. For fo and the harmonic frequencies 
(i.e. 2/o, 3/o, etc.) we observe a multiplet structure of the 
same properties as discussed earlier (see Fig. [3} • 

Even more complicated behaviour is observed for the 
model of T e B = 5876 K. We note that this is the model of the 
longest, peak modulation period (Fig. [5}. In fact the period 
of the modulation cycle, if measured from one epoch of maxi- 
mum pulsation amplitude to the other, may strongly vary for 
this model which is well visible in the bottom panel of Fig. [7] 
For the plotted cycles we estimate it to be 96 d, 95 d, 129 d, 
106 d and 107 d. The modulation is quasi-periodic. The first 
return map for this model is presented in Fig. [6] Analysis 
of the inset shows that the radii maxima either follow the 
route along the figure '8' or follow only its one part. As long 
as we continued the model integration (20 000 fundamen- 
tal mode pulsation cycles, i.e., more than 400 modulation 
cycles) no regularity could be spotted in the described be- 
haviour. We also note that consecutive radii maxima do not 
follow a single curve, which is well visible in the respective 
inset in Fig. |S] Consequently the modulation pattern varies 
from one modulation cycle to the other. This quasi-periodic 
behaviour is also detected in the frequency spectrum, which 
is plotted in the bottom panel of Fig. [5] At fo and its har- 
monic we detect the modulation multiplet with the sepa- 
ration corresponding to f m . The inverse, /m 1 , describes the 
mean modulation period, which is equal to 106.7 d. At the 
sub- harmonic frequencies broad bands of peaks with no clear 
structure are visible. These peaks are not coherent just as are 
the modulation peaks around fo and its harmonics, clearly 
indicating the chaotic nature of the model. Still the modu- 
lation can be described as more or less periodic, as it is the 
case in several Blazhko RR Lyrae stars, with strongly d iffer- 
ing modulation cycles (e.g. iGuggenberger et al.ll201ll ). We 
note that properties of Fourier spectrum for the discussed 
model closely resemble the properties of Fourier spectra for 
Blazhko RR Lyrae stars showing the period doubling ef- 
fect. For these stars bunches of peaks are present at half- 
integer frequencies indic ating the irregular nature of period 
doub ling phenomenon (|Szab6 et al.l |2010| ; iKolenberg et al.l 

i53). 

As the models' temperature is increased beyond 5876 K, 
the figure '8' is broken into two separate loops (Fig. [6]). 
Strictly periodic modulation, as described in Section 12.11 
is observed again. 

Analysis of the radius variation and return maps for our 
models provide a hint why the modulation period varies for 
these models (as depicted in Fig. [5} and why the modula- 
tion may become quasi-periodic. We focus attention on the 
three previously discussed models (of T o ff=5860 K, 5872 K 
and 5876 K) , for which the radius variation is plotted in 
Fig. [7] During the modulation the amplitudes of the excited 
modes, i.e. of the fundamental mode (A(t) in the follow- 
ing), and resonantly coupled first overtone (B(t)) vary. For 
the three models we observe that the period doubling is 
quenched at some phases of the modulation cycle. It is visi- 
ble in the return maps (radii maxima fall close to the diago- 
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Figure 6. First return maps for maximum radii (in units of R©) for several models of 128 Lq. The insets in three panels show the zoom 
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Figure 7. Radius variation for the three models of 128 Lq and different effective temperatures. Right panels show the zoom into a 
100-day segment. 



nal, see insets in Fig. [6]) and clearly in Fig. [7] The amplitude 
of the overtone mode must be very small then. We conclude 
that the model evolves in direct proximity of the point cor- 
responding to unstable single periodic pulsation, i.e. a point 
of A — const and B = 0. This unstable limit cycle corre- 
sponds to a single point on the diagonal of the return map. 



The topology of return map for T e ff =5876 K (Fig. [6} suggest 
that it is a saddle point. The exact location of the unsta- 
ble single periodic limit cycle co uld be determined with the 
help of the relaxation technique (IStellingwerJll974 ). Unfor- 
tunately it is not implemented in our code. Nevertheless, in 
the next Section, with the help of the amplitude equations 
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Figure 8. Frequency spectra after prewhitening with the fun- 
damental mode and its harmonics for two models of 128 Lq se- 
quence. Note different scales on the vertical axis in the left and 
right sides of the Figure. Removed frequencies are marked with 
dashed lines, while arrows (top panel) show the exact location of 
the sub-harmonic frequencies, l/2/o and 3/2/q. 



formalism, we show that indeed the observed properties of 
the modulation result from the dynamical evolution close 
to the unstable saddle point. We will show that at a saddle 
point the trajectories describing the system evolution are at 
first attracted, but repelled later on, and that evolution is 
very slow (the closer to saddle point, the slower the system's 
evolution). Thus, the models evolving closer to the unstable 
saddle point have much longer modulation periods (and a 
much longer phase of suppressed period doubling). In addi- 
tion once the trajectories evolve very close to the saddle, the 
system becomes chaotic. When the system evolves far from 
the saddle point the modulation is strictly periodic and oc- 
curs on a shorter time scale. 

The described properties of radius variation are qual- 
itatively the same for our least luminous model sequence 
(120 Lq). In Fig. |9]we plot the return maps for three models 
located at the peak of modulation periods (Fig. [5}. Again 
we observe that the radii maxima fall close to the diagonal 
along a broad band (see insets for models of 5850 K and 
5851 K). For the model of T a g — 5851 K the radii maxima 
may follow a route along the loops which approach the diag- 
onal (and the presumed saddle point) at different distances. 
Corresponding modulation cycles have different length. 

In Fig. [5]we also note a rapid increase of the modulation 
periods at the cool side of the modulation domain, particu- 
larly pronounced for the L — 120 L0 models. This increase 
again can be explained by the models' evolution close to the 
unstable saddle point. 

Finally we note that also for 136 L.0-sequence the in- 
crease of modulation period is visible (Fig. [5]) which is, how- 
ever, not as prominently marked as in case of the lower lu- 
minosity models. 

We now turn to a discussion of the described modula- 
tion properties using the amplitude equations' formalism. 



3 AMPLITUDE EQUATIONS 

The dynamical behaviour we have found in our mod- 
els can be well captured with the a mplitude equations 
forma lism. We follow the analysis of iBuchler fc Kollath] 
l|201ll ). who considered the 9:2 (9ojo = 2ujg) resonance 
responsible for t he period doubling d e tected recently in 
RR Lyrae stars jKolenberg et all l20iri ISzabo et all l201fj| ; 
iKollath. Molnar fc Szabd 2011 ). All these stars show the 
Blazhko effect, m ore or less periodic modu lation of the pul- 
sation amplitude. IBuchler fe Kollathl |201ll ) showed that the 
same resonance is able to cause the amplitude modulation. 
They conducted the parametric study of the appropriate 
AEs and besides the fixed points they identified two types 
of oscillatory solutions, either strictly periodic or irregular. 
Here we follow their approach to study the 3:2 resonance be- 
tween the fundamental mode and the first overtone which is 
responsible for th e period doubling e ffect in BL Her hydro- 
dynamic models (|Smolec et al.ll2012| ). As we noted in Sec- 
tion 12.11 the stable period two cycle loses its stability and 
gives rise to a stable cycle in which pulsation is modulated, 
either periodically or in a quasi-periodic manner. Below we 
show that this behaviour can be captured with the ampli- 
tude equations. 

The appropriat e set of complex amplit ude equa- 
tions was given by Moskalik fe Buchlerl |l99y) (see also 
iMoskalik fe Buchlerll 19931 ): 

— = (iiu a + K, a + Q a \a\ 2 + T a \b\ 2 )a + c a a* 2 b 2 , (la) 

-3£ = (iu b + K b + Q b \b\ 2 +T b \a\ 2 )b + c b a :i b* . (lb) 

The non-resonant part is truncated at the third order. Qk 
and Tfc are cubic, self- and cross-saturation coefficients, re- 
spectively. The last terms on the r.h.s. of ([TJ are resonant 
coupling terms, with resonant coupling coefficients c a and c b . 
These, together with complex amplitudes, are represented as 
a = Ae i0 °, b = Be i0b , c a = C a e iSa and c b = C b e iS " , which 
allows to write the above AEs as set of three first order 
differential equations with real coefficients: 

dA 



= (Ka+^(Qa)A 2 +U{T a )B 2 )A + 
+ CaA 2 B 2 COs(T + S a ), 

= (K b + M(Q b )B 2 + $t{T b )A 2 )B + 

+C b A :i B cos(r-<5 6 ), 

= 2A+(2Q(T b )-3Q(Q a ))A 2 + 

+ (2Z(Q b )-3Q{T a ))B 2 + 

-ZCaAB 2 sin(r + S a ) - 2C b A 3 sin(r - 8 b ) . 

Above F is a combination of mode phases, F = 2<j> b 
and A is a frequency mismatch parameter, 



dt 

dB 
~dt 

dT 
~dt 



A 



: LO b 



3w a /2. 



(2a) 
(2b) 
(2c) 
)• 

3</>a, 

(3) 



Exact analytical solutions of amplitud e equations are pos- 
sible o nly for strongly simplified cases. IMoskalik fe Buchlerl 
( 1990) considered the AEs for half-integer resonances be- 
tween the fundamental mode and the higher-order over- 
tones, (2n + 1)ojo = 2uj k , neglecting the cross-saturation 
and resonant coupling term in the equation for the domi- 
nating fundamental mode (i.e., T a — and c a = 0), as well 
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as neglecting the imaginary parts of the cubic coupling co- 
efficients. Such system can have two non-trivial fixed points 
(dA/dt = dB/dt — dT/dt — 0), either single- mode fixed 
point with non-zero amplitude of the fundamental mode 
only (A = \J— K a /^t{Qa) and B — 0) or two-mode fixed 
point with non- vanishing amplitudes of both modes (that in 
general have to be derived numerically). Because of the half- 
integer resonance, such solution corresponds to a period- 
doubled limit cycle. Oscillatory solutions for the amplitudes 
are possible as well and were found during the numerical 
integr ation of eqs. ([2]). To sim plify the equations we fol- 
lowed iBuchler fe Kollathl (|201lh and neglected the imagi- 
nary parts of the cubic saturation coefficients. Periodic mod- 
ulation of amplitudes, corresponding to hydrodynamic mod- 
els described in Section [2] was relatively easy to find. We 
adopted the following set of the saturation and coupling co- 
efficients: K a = 1.0, K b = -0.01, 5R(Q a ) = 1, K(Q 6 ) = 1.1, 
5R(T a ) = 5.0, U(T b ) = 1.0, C a = 5.75, C b = 4.0, 8 a = 1.95, 
5 b = 1.0 and integrated the AEs with different mismatch pa- 
rameter A. These models correspond to our hydrodynamic 
models of different effective temperatu res. The radius vari- 
ation was reconstructed according to ( Moskalik fe Buchlerl 
1 19901 ) 

SR{t) = Ae ioJat + £ e i(r/2) e i(3/2)c a t + h Q t ^ ^ 

with higher order terms (h.o.t.) neglected. For the pulsa- 
tion frequency, ui a , we choose ui a = 2n/P a = 27r/0.1 so that 
the relative growth rate, K a P a — 0.1 is typical for BL Her 
models. We note, that we do not expect a 1:1 equivalence 
between the hydrodynamic models and results of AEs inte- 
gration, not only because of arbitrarily assumed coefficients 
of eq. but also because of neglect of the higher order 
terms in radii reconstruction (eq. 2]). 

In Fig. [10] we show the radius variation for three dif- 
ferent values of A. Qualitatively the same behaviour is vis- 
ible as in the hydrodynamic models (Fig. [7} - clear ampli- 
tude modulation on top of the period doubling behaviour. 
The modulation is periodic for A = 0.4 and A = 0.8, and 



is clearly quasi-periodic for A = 0.fJ3- The corresponding 
return maps on successive radii maxima are presented in 
Fig. [TTJ Diamond and triangles mark the location of radii 
maxima corresponding to the single-mode fixed point (di- 
amond) and two- mode fixed point (triangles). The fixed 
points are unstable as their stability eigenvalues indicate. 
For the single-mode point the three eigenvalues are all real, 
one is positive and two are negative, thus the single-mode 
fixed point is an unstable saddle point. For two-mode fixed 
point one eigenvalue is real and negative, the other two are 
a complex conjugate pair with positive real parts. The two- 
mode fixed point is thus an unstable spiral point. For mod- 
els with A = 0.40 and A = 0.80 the radii maxima fall along 
single curves encircling the unstable two mode fixed points 
and far from the unstable single- mode point. These return 
maps are qualitatively the same as for our hydromodels that 
display the strictly periodic modulation (e.g. the models of 
T cff = 5844 K or T cff = 5872 K in Fig. |BJ. The behaviour 
is very different for A = 0.6. The radii maxima may fall at 
direct proximity of the unstable single-mode point and may 
follow a very different routes during the modulation. This 
is qualitatively the same behaviour as described for models 
of r e a = 5876 K (128 L ) and of T e g = 5851 K (120 Lq) in 
Section ^. 21 It results from the evolution at the proximity of 
the unstable saddle point. 

To analyse the dynamic s of ou r models in more detail 
we follow [Buchler fc Kollathl |201lf ) and consider the phase 
plots of amplitudes, i.e. plots of A(t) vs. B(t), which are 
displayed in Fig. for the discussed cases. Note that the 
(A, _B)-plane is a subset of the full 3D phase-space (A, B, 
r) and that now we focus on the long term variation, with 
pulsations filtered out. For A = 0.4 and A = 0.8 the phase 
plots are very simple - the trajectory is a single loop and 
resulting modulation is strictly periodic. In the middle panel 
of Fig. [T5] we plot the trajectories for the most interesting 
case of A = 0.6. Clearly, the trajectory is not a single curve 
now, which is best visible at the proximity of the unstable 

1 The apparently large values of the mismatch parameter, A, 
result from the adopted normalization (n a = 1). 
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Figure 10. Radius variation for three integrations of AEs with different A parameter. Note the different lengths of consecutive modulation 
cycles for A = 0.60. 
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Figure 11. First return maps for radii maxima for models with different A. Diamond and triangles mark the location corresponding to 
unstable single-mode and two-mode fixed points, respectively. 



single-mode point (A = 1, B = 0). The phase plot indi- 
cates the presence of the strange attractor, which may be 
further confirmed through construction of the first return 
map on the consecutive maxima of amplitude A(t). It is 
plotted in Fig. [T3] a nd is very similar to that of t he Lorentz 
attractor (see e.g. iBerge. Pomeau fc Vidallll984 ). albeit, at 
close inspection, it shows the double-cusp feature (see in- 
set) . The slope of the vs - the A^ ax at the intersection 
with diagonal is greater than unity (| slope | > 1) and there- 
fore the behaviour is chaotic. T he character of the return 
map is also very similar to that of iBuchler fc Kollathl (|201ll ) 
(see their fig. 3) who studied the 9:2 resonance suspected 
to cause the Blazhko effect in RR Lyrae stars. Clearly the 
two half-integer resonances lead to similar types of dynam- 
ical behaviour, either to a strictly periodic modulation, or 
to a chaotic one. We stress however that our models, both 
hydrodynamic and based on the integration of AEs are not 
'too chaotic' - the modulation may be described as quasi- 
or more or less periodic. This is exactly what we observe in 
Blazhko RR Lyrae stars. 

The evolution close to the saddle point, in particular 



the differences in periods of the consecutive modulation cy- 
cles, can be understood through analysis of the flow field, 
i.e. (A, B), which varies in time. It is tangent to the trajec- 
tory at the point {Ait), B(t)). The animation attached as a 
supplementary on-line material shows the evolution for two 
trajectories: one located closer to the saddle than the other. 
The snapshot from this animation is plotted in Fig. 1141 It 
is clear that at the saddle trajectories are at first attracted, 
but repelled later on. Also the evolution at the saddle is 
slow (short vectors in the flow field) - the closer to the sad- 
dle the slower the evolution. Thus, the modulation cycles 
corresponding to the two trajectories have different length. 

In Fig. [15] we plot the mean modulation periods and 
modulation amplitudes, as defined in Section f2.ll for a set 
of models with different A. Just as in the case of the hy- 
drodynamic models (of 120 L© and 128 Lq, see Fig. [5]) we 
observe the increase of the mean modulation period, accom- 
panied with the increase of the modulation amplitude. The 
longest modulation period is found for A = 0.6 for which 
the modulation is chaotic. 
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Figure 12. Trajectories in the (A(t), B(t))-plane for several modulation cycles and different A. Unstable single-mode and two-mode 
fixed points are marked with diamond and triangles, respectively. 
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Figure 13. First return map on the successive amplitude maxima 
of amplitude A(t); A^lx vs. A^ ax . 




Figure 14. Snapshot from the animation attached as a supple- 
mentary online material. Trajectories for two modulation cycles 
passing either close or farther away the saddle are plotted. The 
flow field is visualised at the instant marked with the circle. 
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Figure 15. Mean modulation periods and modulation ampli- 
tudes as a function of A. 



4 DISCUSSION 

• Origin of the period doubling and of modulation of pul- 
sation. In Fig. \T\ loci of several half-integer resonances are 
plotted. In principle each of them may cause the period dou- 
bling effect, provided that the resonant mode is not damped 
too strongly, a nd that respective model f alls close to the res- 
onant centre (|Moskalik fc Buchler|[l990h . These conditions 
rule out all the resonances displayed in Fig. [1] except the 
3:2 resonance between the fundamental mode and the first 
overtone. In Fig. [TS]we plot the amplitude of the period dou- 
bling (amplitude of the highest peak at around /o /2) vs. the 
mismatch parameter for the three half-integer resonances. 
Only models of 136 Lq are plotted as for this luminosity we 
conducted a model survey through the full period doubling 
domain. It is clear that only for the 3:2 resonance all the 
models are close to the resonance centre, with | A.a : 2 1 < 0.06 
(where A3.2 = oj\/ujq — 1.5). Considering e.g. the 7:2 reso- 
nance with the fourth overtone (middle panel of Fig. I16p . 
although some models are located at direct proximity of the 
resonance centre, other models within the period doubling 
domain may be located as far as A7 : 2| > 0.5, which rules out 
this resonance as a possible cause of period doubling. The 
same holds for the 9:2 resonance with the seventh overtone. 

The domain with the modulation of pulsation lays within 
the period doubling domain. Its location with respect to the 
mismatch parameters is also plotted in Fig. 1161 in which 
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we used the amplitude of the highest modulation peak in 
the vicinity of fo (fo + fm or fa — / m ) to characterise the 
strength of modulation. For clarity, only the domain for 
L — 136 Lq is plotted, the domains for other luminosities 
fall roughly at the same location. Already from Fig. [1] it is 
clear that the domain with modulation of pulsation follows 
(at lower luminosities) the loci of the 7:2 resonance with 
the fourth overtone. Indeed, in Fig. [16] the resonance cen- 
tre falls in the middle of the modulation domain and all 
of the models lay in close proximity to the resonance cen- 
tre - as one may read from Fig. [16] (|A7 : 2| < 0.05 for all 
the models). This is however true also for the 3:2 resonance 
that causes the period doubling. Although the models are 
at offset with respect to the resonance centre, they all fall 
within 0.02 < A3. 2 < 0.04. As the modulation domain arises 
within the period doubling domain, which is caused by the 
3:2 resonance, it is therefore natural to assume that the same 
resonance causes the modulation. It is confirmed with the 
analysis of t he amplitude equations p resented in Section [3] 
(and also in iBuchler fc Kollatbl l201ll ) . No additional reso- 
nance is needed to cause the modulation. It appears within 
the period doubling domain, not at its centre but at positive 
A fFig. 1151 note a different definition of A at this plot). Note 
that the positive mismatch as defined in eq. in Section [3] 
corresponds to A3 : 2 > 0, just as we have for our hydromod- 
els (Fig. 1 16 P . Taking into account that all the behaviours we 
have found in our hydrodynamic models (period doubling, 
periodic and quasi-periodic modulation) can be reproduced 
with AEs for the 3:2 resonance only, we conclude that this 
resonance is the only cause of both period doubling and am- 
plitude modulation in our hydrodynamic models. 

• Reliability of the computed models. Our models have 
strongly decreased eddy viscosity as compared to the stan- 
dard models of BL Her s tars, i.e. a m = 0.05 to be compared 
with a m = 0.2 we used in ISmolec et ahll|2012h . Consequently 
the pulsation is violent and spurious spikes appear in the 
luminosity curve. Still the models may represent the phe- 
nomena that may occur in real BL Her stars. Certainly it 
is so for the period doubling discovered in these stars only 
recently, 20 years after the effect was predicted based on 
hydrodynamic models (see Introduction). No modulation of 
pulsation was observed in BL Her stars so far. 

The turbulent convection model that we use (|Kuhfufil 
1986) is a very simple ID formula with several free pa- 
rameters. The treatment of viscous effects seem to be the 
weakest point of this and similar models, which is not sur- 
prising. The viscous dissipation of turbulent energy takes 
place over many length scales and its one parameter de- 
scription is, out of necessity, provisional. We note that 
the current convective models cannot reproduce a double- 
mode pulsation unless unphysical neglect of negative buoy- 
ancy below the convective zone leads to artific ial viscous 
damping there jSmolec fc Moskalik|[2008b1 . l2010t ). With the 
purely radiative models double-mode pulsation couldn't be 
reproduced either, unless the artificial dissipation (replaced 
in convective codes by v iscous dissipation) was decreased 
(|Kovacs fc Buchlerlll993 | ). Similarly, the period doubli ng in 
RR Lyrae models of iKollath. Moinar fc Szabd |201lf ) was 
produced assuming lower values of eddy viscosity. Thus, dis- 
covery of modulation of pulsation in BL Her stars in the 
future cannot be excluded. 

• Models with chaotic modulation of pulsation. In case of 
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L — 136 Lq we computed additional models across the in- 
stability strip. This allowed us to detect another domain 
in which pulsation is modulated, this time however, in a 
very chaotic manner. In these models we observe a wealth 
of dynamical behaviours. In particular we identify several 
stability windows within the chaotic regime, with stable pe- 
riod k cycles (e.g. with k — 3, k — 7 or k = 9). They arise 
due to the tangent bifurcation and, as the control parameter 
(T e ff) is increased, undergo a series of pitchfork bifurcations 
(sub- harmonic cascade). These models will be subject of a 
separate paper (Smolec & Moskalik, in prep.). 
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5 IMPLICATIONS FOR THE BLAZHKO 
EFFECT IN RR LYRAE STARS 

The modulation of pulsation is a common property of 
RR Lyrae stars, where it is called the Blazhko effect. 
In several stars period doubling was also detected and 
was reproduced by hydrodyna mic RR Lyrae models of 
iKollath. Molnar fc S zabo ( 2oTl|) , which however, show no 
trace of modulation. Period doubling in the hydrodynamic 
RR Lyrae models is caused by the 9:2 resonance between 
the fundamental mode and the ninth overtone . The a nalysis 
of amplitude equations bv lBuchler fc Kollathl l|201ll 'l shows 
that the same resonance may cause the modulation of pul- 
sation. Confirmation of this result through direct hydrody- 
namic modelling would strongly support this radial reso- 
nance model. Our hydrodynamic BL Her models are the 
first models that demonstrate that the mechanism may in- 
deed work. Except for larger luminosities, BL Her stars are 
siblings of RR Lyrae stars. The resonance that causes the 
modulation of pulsation in our models, 3:2 resonance with 
the first overtone, cannot work in case of RR Lyrae stars. 
Its locus falls beyond the RR Lyrae instability strip. It is 
clear however, that the underlying dynamics is the same for 
the 3:2 and for the 9:2 resonanc es as studies of amplitud e 
equations indicate (Section l3l and IBuchler fc Kollatbll201l[ ). 

The great advantage of the radial resonance model 
is that it can cause the irregular modulation, observed 
in Blazhko RR Lyrae stars. The striking evidence comes 
from the near ly continuous observations of the satellite mis - 
sions (see e.g. lBenko et alj|20ld : ICuggenberger et al.ll2012l) . 
Our models together with analysis of amplitude equations 
demonstrate how these irregularities may arise in the frame- 
work of the resonance model. For some model parameters 
the modulation is no longer strictly periodic, but becomes 
chaotic. The models indicate that it happens when ampli- 
tudes of the involved fundamental and first overtone modes 
approach the unstable saddle point which corresponds to 
single periodic fundamental mode pulsation. Then the pe- 
riod doubling is quenched and the variation of the modal 
amplitudes is very slow: the closer to the saddle, the slower 
the modulation. Thus consecutive modulation cycles dif- 
fer. It is worth noting that our models showing the ir- 
regular modulation of pulsation reproduce other features 
of Blazhko RR Lyrae stars that simultaneously show the 
period doubling effect. We observe qualitatively the same 
properties of Fourier spectra at half integer frequencies - 
bunches of peaks, indicating the irregular nature of period 
doubling phenomenon, as well as swapi ng of higher/lower 
maximum cycles duri ng the modulation |Szabo et al.ll20l"ol : 
iKolenberg et al.ll201ll . Molnar et al., submitted). 

Comparison of the properties of modulation in our 
BL Her models to that observed in the Blazhko RR Lyrae 
stars points to the possible challenges for the radial res- 
onance model. One of them is a very strong period dou- 
bling effect in all our models. On the contrary, only in a few 
Blazhko stars observed with satellite mis sions the effect wa s 
detected and its amplitude is rather low l|Szabo et al.| [2010L 
The strong period doubling is also evident in our models 
based on the amplitude equations integration, as well as 
in models corresponding to the 9:2 resonance analysed by 
IBuchler fc Kollathl l|201lf ). which should be operational in 
RR Lyrae stars. The lack of period doubling in most of the 



Blazhko stars indicate that the modulation affects predom- 
inantly the fundamental mode, while the amplitude of sus- 
pected resonant mode must remain quenched during most 
of the modulation cycle. Our models indicate that the pe- 
riod doubling is quenched while the model evolves close to 
the single-mode saddle point, where the amplitude of the 
resonant mode approaches to zero. The resonant mode will 
always be modulated in the resonant scenario - the modal 
amplitudes must encircle the unstable two mode fixed point 
- see Fig. 1121 To decrease the modulation of the resonant 
mode and hence to quench the period doubling, one has to 
bring the two-mode fixed point closer to the B(t) = axis. 
An extensive parameter study of the amplitude equations is 
needed to point under which conditions such modulation is 
possible. 

Irregularities are often detected in the Blazhko stars. 
In our models they seem to be restricted to a very nar- 
row parameter range, most of our models show a strictly 
periodic modulation. This might not be a problem for the 
resonant scenario if one takes into account the dispersion 
of stellar parameters (masses, luminosities, metallicities) of 
real stars or if one allows other resonances to play a dy- 
namical role. The latter possibility is worth considering. 
Nearly 50 per cent of RR Lyrae stars (pul sating in the 
fundamental mode) show the Blazhko effect jjurcsik et al.l 
2009; iBenko et al.l |2010T ). It seems unlikely that the same, 
high-order resonance operates in all these stars. All other 
resonances within the radial modes face the same prob- 
lems (or even harder, as most of the overtones are heavilly 
damped). The half-integer resonances with non-radial modes 
may present a solution. Non-radial modes are densely packed 
in between the radial modes, thus there is always a candi- 
date for resonant coupling. On the other hand, the mode 
selection would be a crucial problem here. We note that 
in recent years non-radial modes were detected in several 
RR Lyrae star s, both in the ground-based and space ob- 
servations (e.g. Gruberbauer et al1l2007l ; lOlech fc Mo skalik 
120091; ICuggenberger et alj|2012h . 
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SUPPORTING INFORMATION 

Additional Supporting Information may be found in the on- 
line version of this article: 

Animation. Model's evolution in the B(t) vs. A(t) plane 
(A = 0.6) for the two modulation cycles. Diamond marks 
the location of unstable single- mode fixed point (saddle). 
Circle shows the model's position during the evolution. The 
instantaneous flow field (A(t), B(t)) is plotted with vectors. 
Snapshot from the animation is plotted in Fig. 1141 



